Phase Transitions in Pressurised Semiflexible Polymer Rings 
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We propose and study a model for the equilibrium statistical mechanics of a pressurised semiflex- 
ible polymer ring. The Hamiltonian has a term which couples to the algebraic area of the ring and 
a term which accounts for bending (semiflexibility) . The model allows for self-intersections. Using a 
combination of Monte Carlo simulations, Flory-type scaling theory, mean-field approximations and 
lattice enumeration techniques, we obtain a phase diagram in which collapsed and inflated phases 
are separated by a continuous transition. The scaling properties of the averaged area as a function 
of the number of units of the ring are derived. For large pressures, the asymptotic behaviour of the 
area is calculated for both discrete and lattice versions of the model. For small pressures, the area 
is obtained through a mapping onto the quantum mechanical problem of an electron moving in a 
magnetic field. The simulation data agree well with the analytic and mean-field results. 

PACS numbers: 64.60.Cn,05.70.Fh,05.50.+q,05.10.Ln 



INTRODUCTION 



Fluid vesicles obtained via the self-assembly of am- 
phiphilic molecules exhibit a variety of shapes in thermal 
equilibrium. Such shapes can be understood in terms 
of the energy minimising configurations of a curvature 
Hamiltonian, under the constraints of fixed enclosed vol- 
ume and surface area [11,11,11. Shape Chan ges arise when 
solutions of the Euler-Lagrange equations representing 
distinct shapes exchange stability. However, the non- 
linearity of these equations, if no special symmetries are 
assumed, necessitates purely numerical approaches. Fur- 
ther, while the curvature modulus in bilayer lipid mem- 
brane systems is often large, so that thermal fluctua- 
tions about the minimum free energy structure may be 
ignored, the more general problem of understanding the 
thermodynamics of such shape transitions is a formidable 
oneQ. 

The two-dimensional version of the vesicle problem 
is a polymer ring of fixed contour length, whose en- 
closed area A is constrained through a coupling to a 
pressure difference term p. Leibler, Singh and Fisher 
(LSF) [5] performed a Monte Carlo and scaling study of 
two-dimensional vesicles, modelled as closed, planar, self- 
avoiding tethered chains, accounting for both pressure 
and bending rigidity. In this model, the ring polymer is 
obtained by connecting the centres of impenetrable par- 
ticles of fixed radius with tethers of a fixed maximum 
length, while enforcing self-avoidance. LSF showed the 
existence of a phase transition at p = 0, separating a 
branched polymer phase for p < from an inflated phase 
for p > 0. At the transition point, the ring is described by 
a self avoiding polygon. Various fractal and non-fractal 
shapes that arise in these models have also been investi- 
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Analytic studies of this class of models present many 
difficulties, arising principally from the self-avoidance 
constraint. Nevertheless, the relatively simple structure 
of the LSF model has stimulated a considerable body of 
work, largely in exact enumeration studies of lattice ver- 
sions of the original continuum model and its variants 
[i, lal lid, llll, IJJ, 1131, IH 1151, IJJI- Most of these studies 
have concentrated on the behaviour of the system in the 
thermodynamic limit in the region p < Q. However, the 
p > case can exhibit interesting crossover behaviour for 
large but finite systems. 

The consequences of relaxing the self-avoidance con- 
straint were studied in Refs. [ij lH, H^. In the models 
studied in these papers, the ring was allowed to inter- 
sect itself, with the pressure term coupled to the alge- 
braic area [13, HB] or to its square 19]. The particles 
linked to form the polymer were coupled through har- 
monic springs [13, 13 , thus allowing for the extensibility 
of the chain. We shall refer to this model as the Extensi- 
ble Self-Intersecting Ring (ESIR). The ESIR model can 
be solved exactly. The solution yields collapsed and in- 
flated phases of the ring separated by a continuous phase 
transition that occurs at a critical value of an appro- 
priately scaled pressure ^Ts'l. However, the model has a 
major shortcoming in that the inflated phase is an un- 
physical one in which the ring expands to an infinite size. 
In a more realistic model, such an expansion would be 
limited by the finite size of individual link lengths. 

The unphysical nature of the inflated phase in the 
ESIR model has been addressed in recent work[23|, in 
which particles are joined by bonds of flxed length, as 
opposed to springs. The Hamiltonian has a term where 
the pressure couples to the algebraic area, as in the 
ESIR model. The transition survives as a continuous 
phase transition with mean-field exponents, separating 
collapsed and inflated regimes of the ring. We shall refer 
to this model as the Inextensible Self-Intersecting Ring 
(ISIR). 

The model proposed in this paper incorporates a bend- 



2 



(a) 




(c) 


(d) 

C ) 



FIG. 1: The collapsed to inflated phase transition as the 
pressure is increased. The different panels correspond to (a) 
J = 0, p < Pc; (b) J = 2, p < pc, (c) J = 0, p = pc, and (d) 
J = 0, p > Pc- 



ing energy into the ISIR model along standard lines for 
semiflexible polymers. We retain the coupling of the 
signed pressure to the algebraic area noting, as argued 
in [l^l, that this difference, while vastly increasing the 
tractability of the problem, makes little difference to com- 
putations within the inflated phase. 

The continuum problem we address is the following. If 
the polymer chain is specified by the curve r(s), where s 
is the arc-length along the curve, we consider the Hamil- 
tonian 



^ = —2 



ds 



df\ K 
r X — ] ■ z + - 

as 



ds 



ds^ 



(1) 



Here, ds{f x dr/ds) ■ z is the infinitesimal (signed) area 
which must be integrated over the internal variable s to 
obtain the total algebraic area. The quantity L — Na is 
the length of the polymer ring of N units where a is the 
size of the basic monomer. The continuum limit is taken 
such that N ^ oo and a ^ 0, keeping L = Na fixed. 
The parameter p = pin — Pout represents the pressure 
differential between the inside and the outside of the ring 
and K is the bending rigidity of the chain. We measure 
energies in units of k^T. The inextensibility constraint 
is imposed through 



= 1, 



(2) 



where t{s) is the unit tangent vector. When p — 0, this 
model reduces to the worm-like chain model, with con- 
figurations constrained by the closure requirement. 



We will work with two discretized versions of the above 
Hamiltonian. The first has N particles in the continuum 
connected through fixed length links, forming a ring poly- 
mer whose equilibrium configurations are constrained by 
the pressurisation and bending energy terms above, but 
where the ring can intersect itself at no energy cost. We 
will refer to this version of our model as the "discrete 
model". The second is a square lattice version of the 
same problem with the particles constrained to lie on the 
vertices of the lattice. We will refer to this version as the 
"lattice model" . We discuss the differences and similari- 
ties between the two versions. 

We use a combination of analytic and numerical meth- 
ods to study these models: Flory type scaling theory 
for the scaling of the area as a function of pressure, 
Monte Carlo simulations for different pressures and bend- 
ing rigidity, mean field approaches and exact enumera- 
tion. 

In Fig. [1] we show typical configurations obtained from 
Monte Carlo simulations of the discrete model in four 
limits. These are configuration snapshots across the col- 
lapsed to infiated phase transition, for different values of 
the bending rigidity J of the discrete model (J oc k), 
as the pressure p is varied. Fig. [Ha) shows the collapsed 
phase for the case where the bending energy is zero, while 
Fig. [Ijb) illustrates a typical ring configuration at an in- 
termediate value of the bending rigidity, but still within 
the collapsed regime. In Fig. [IJc), we show a typical 
configuration close to the transition between collapsed 
and inflated phases. Last, Fig. [Ud) illustrates the fully 
inflated ring. 

We summarise our main results below. We show that 
there is a continuous phase transition in the scaled pres- 
sure p {Np/Ai:) - bending rigidity (J) phase diagram, 
which separates a collapsed phase in which area oc N ^ 
from an inflated phase in which area oc iV^ (see Fig. [2|) . 
The phase boundary for the discrete model is obtained 
as Pc = [/o(J) -/i(J)]/[/o(J) +/i(J)], where the /(J)'s 
are modified Bessel functions. For the lattice model, the 
phase boundary is obtained as Pc = e~'^ . 

These results are obtained by first solving the J — 
case exactly and then incorporating the effects of a 
nonzero J through a scaling argument. For the collapsed 
phase, the free energy for nonzero J is calculated by the 
same method. In the inflated regime, we resort to mean 
field theories. We employ two types of mean-fleld theo- 
ries: In the flrst, the inextensibility constraint is satisfied 
exactly but the closure condition is satisfied only on av- 
erage. In the second, we impose the closure condition 
exactly but satisfy the inextensibility constraint only on 
average. The dependence of the area on p for p oo 
is calculated. The behaviour near the transition line is 
obtained through a Flory type scaling theory. 

The rest of the paper is organised as follows. In Sec. |TT] 
we define our models more precisely. Section [ITT] con- 
tains the details of the numerical methods used, includ- 
ing the Monte Carlo and exact enumeration algorithms. 
In Sec. IIVI we discuss a Flory-type scaling theory valid 
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FIG. 2: The phase boundary between collapsed and inflated 
phases for a semi-flexible polymer ring as obtained by two 
different methods, a scaling analysis based on Flory-type ar- 
guments and mean-field theory. 

for the semi-flexible case. Section |V] describes mean- field 
approaches to this problem: (a) a simple density-matrix 
based single-site mean-field approach which captures the 
properties of the inflated phase to very high accuracy 
but is inadequate for the coUapsed phase and (b), a less 
accurate harmonic spring mean-field theory, which is ca- 
pable of describing both collapsed and inflated phases, 
fn Sec. IVIl we discuss the behaviour around the critical 
point in greater detail. Our Sec. IVIl] contains results for 
the asymptotic behaviour of the area as well as a descrip- 
tion of the appropriate scaling function for the area in the 
lattice case, as a function of N. Section IVlIII contains a 
summary and conclusions. In Appendix \X\ we present a 
solution of the problem when J = by drawing an anal- 
ogy with the quantum mechanical problem of the motion 
of an electron in a magnetic field. 



lowing standard procedures as 

N 

Hb = -JY,U-U+1, (5) 

i=l 

where the bending rigidity J of the discrete model is pro- 
portional to the continuum bending rigidity and t is the 
unit vector in the direction of t. The inextensibility con- 
dition is imposed through 

\n - f,-i\ = \t,\ ^ a ^ I. (6) 

Since the tangent vectors have unit norm, we can repre- 
sent them as ti = {cos 6i, sin 9i), where 9 e [0, 27r). In 
terms of these variables, the partition function is 

i j=0 \k=Q ) 

(7) 

On a square lattice, the model remains essentially the 
same except for restrictions on the angles 6i. Now Qi is 
allowed to take values 0, 7r/2, tt, 37r/2 such that all the 
particles are on the vertices of the square lattice. 

III. NUMERICAL METHOD 

In this section, we describe the numerical methods 
used. For the discrete version of our model, we use Monte 
Carlo simulations (described in Sec. IIII Ap while for the 
lattice problem on the square lattice, we use an exact enu- 
meration scheme (described in Sec. IIIIBj) . The analytic 
results we obtain for our model, described in later sec- 
tions, provide useful benchmarks for the numerical work. 



II. MODEL 

Consider a closed chain of N monomers in two dimen- 
sions. Let the positions of the particle be denoted 
by the vector and the corresponding tangent vectors 
by tj — Tj+i — Tj, 2 = 1,2,..., TV. For a closed ring, 
rjv+i = n., or equivalently, X)i ^» ^ 0- The algebraic or 
signed area area enclosed by the ring is given by 

^ N N J-l 

i=l j=l k=l 

As can be either positive or negative. 

Coupling this algebraic area to pressure, we obtain the 
energy term. 

Hp = -pAs. (4) 
The bending energy cost can then be written down fol- 



A. Monte Carlo Simulations 

The algorithm for the Monte Carlo simulation of the 
discrete model consists of two basic moves [2ll[22j: a sin- 
gle particle flip and a global flip. In the single particle 
flip, a particle is picked at random and reflected about 
the straight line joining its two neighbours (see Fig.[3fa)). 
The move is accepted using the standard Metropolis al- 
gorithm. Since the energy computation involves only 
nearby sites, the move is efficient and fast. In the global 
flip, two particles of the ring are chosen at random and 
the section of the ring between them is reflected about 
the line joining the two particles (see Fig.[3]Jb)). The en- 
ergy calculation now involves 0{N) particles and is thus 
computationally expensive. However, the global move is 
crucial to the study of the case where J 0, since single 
particle moves alone are insufflcient for equilibration in 
this case. 

In the simulations, one Monte Carlo step is defined as 
one global move and N single particle moves made by 
selecting at random particles to be updated. This step is 
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FIG. 3: A schematic representation of the Monte Carlo moves: 
(a) single flip; and (b) global flip. 

then repeated until the system equihbrates. Thermody- 
namic quantities are measured from averages taken over 
independent configurations in equilibrium. 

The initial configuration was chosen to be a regular 
N-sided polygon, but we verified that random configura- 
tions also gave the same results. We performed Monte 
Carlo simulations across a range of pressures for different 
values of J and system size. The system size varied from 
A'' = 64 to TV = 2000. Typically each parameter value 
was run for 4 x 10^ Monte Carlo steps. We waited typi- 
cally for 10^ steps for equilibration, averaging data over 
the remaining steps using independent configurations. 

I 



Similar recursion relations will hold for T/v_|_i(x, ?/; a; -I- 
TN+i{x,y;x,y - 1) and TN+i{x,y;x,y + 1). 
The partition function for the polymer problem can be 
expressed as a sum over areas and bends consistent with 
a given value of the area, i.e., 

2Ar = T^(0,0) = ^CAr(AS)eP^+-^^, (11) 
A,B 

where Cn{A,B) counts the number of closed paths of 
area A in a walk of length N which have B bends. 

We count up to = 150 for different values of J. 
The only limiting factor in going to larger N values is 
computer memory. 



B. Exact enumeration 

We first describe the algorithm for the case J — 0. 
Consider a random walk starting from the origin and tak- 
ing steps in one of the four possible directions. For each 
step in the positive (negative) a;-direction, we assign a 
weight e~^y where y is the ordinate of the walker. 

Multiplying these weights, it is easy to check that the 
weight is e^^ for a closed walk enclosing an area A. 

Let Tn{x, y) be the weighted sum of all A^-step walks 
from (0, 0) to {x, y). It then obeys the recursion relation, 

TN+i{x,y) - e-^^Tw(a:-l,y) + e^^rjv(a: + l,?/) 

+TN{x,y-l)+TN{x,y + l), (8) 

with the initial condition 

7o(a;,?/) = 4,o^y,o- (9) 

Finally, T/v(0,0) gives the partition fimction of the ring 
polymer on a lattice. 

For the semiflexible case, the recursion relation given 
above must be modified, since the ring is no longer a 
simple random walk but a walk with a one step mem- 
ory. We convert it into a Markov process as follows. Let 
Tn{x, y; x' , y') be the sum of weights of all walks reach- 
ing {x,y) in N steps but having been at {x',y') at the 
previous step. These T/v's are now a Markov process and 
depend only on T^-i^s. The recursion relations are then 
straightforward to write down. Rather than give all the 
recursion relations, we provide a representative example 



(10) 

I 

IV. FLORY-TYPE SCALING ANALYSIS 

Flory type scaling theory provides a useful tool to cap- 
ture the scaling behaviour of systems whose free energy 
reflects a competition between two or more terms. Such 
a sc aling theory was proposed for the ISIR model in 
Ref. [2^. A transition from a collapsed to an inflated 
state was predicted to occur at a critical value of the pres- 
sure, whose magnitude scaled with system size as N~^. 
We show how these arguments may be extended to the 
semiflexible case, deriving expressions for the change in 
the critical point and scaling as a function of the bending 
rigidity. 

The free energy consists of three terms describing (i) 



TN+i{x,y;x - l,y) = e ^^[TNix - l,y;x - 2,y) + e^''TN{x - l,y;x,y) 

+e'-'TNix - l,y,x - l,y + 1) + e-^TN{x - l,y;x - l,y - 1)]. 
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FIG. 4: Area versus pressure curves for three J values for 
p < Pc- The points correspond to the discrete case while 
the lines correspond to the lattice case. The inset shows the 
collapse when the curves are scaled as in Eq. (|14|l . The f{x) 
curve in the inset represents the scaling function of Eq. (|16[). 
The data is for iV = 100. 



the entropy of the ring, (ii) the pressure differential and 
(hi) inextensibility of the bonds. When J = 0, these 
terms were argued to be R'^/N, -PR^ and R'^/{'iN^) for 
a ring of size R [20j|, where for the second term it was 
assumed that the area {A) scales as R^ . With semiflex- 
ibihty, we show that a similar scahng form holds except 
for J dependent prefactors. Thus, the free energy takes 
the form 



F 



entropic 



pressure 



F 



inextensibility i 



AttR^ . , , /3( J)i?4 

[a{J)-p\ + " 



N 



Ar3 



(12) 



where we have defined p = Np/An, and a and (3 depend 
on J. 

It is easily seen that a system described by such a Flory 
theory undergoes a continuous transition when the sign 
of the R? /N term changes sign. This occurs at a critical 
scaled pressure Pc{J) which varies with J as 



Pc{J) _ a{J) 
pM " a(0) ■ 



(13) 



When p < Pc{J), then the area follows random walk 
statistics with (A) ^ N. In this regime the R'^/N^ term 
is not important. For length scales much larger than 
the persistence length, the problem is effectively one of a 
freely jointed ring, but with a suitably defined N. Thus, 
we conclude that 



{A{J,N,p)) ^ 



■/ 



P 



N 



Pc{J) \Pc{J) 



P<Pc 



(14) 



where f{x) is a scaling function. The scaling function 
f{x) and Pc can be determined by solving the J — case 
(see Appendix [A| . This gives 



Pc = 47ra(J), 



(15) 



FIG. 5: Comparison of the area ratio {A{J)) / {A{0)) at the 
critical point with the scaling prediction (see Eqns. (|17|l ') for 
the lattice (Eq. (|32])) and discrete (Eq. ((27l)) models. The 
scaling prediction is satisfactory for small J but deviates away 
as J increases. 



and 



J., . 1 cot(7rx) 
fix) ~ 



(16) 



Numerical confirmation of Eqs. (|14p and (fTB|) is provided 
in Fig. m The inset shows that the curves for different J 
collapse onto a single curve when scaled as in Eq. (|14p . 

When p = Pc, the scaling is determined by the R'^/N^ 
term. Thus, {A) - N^/^ / y/j3{J). Thus, 



{A{J)) 

{Am 



i m 



(17) 



To test this relation, we compare the Flory prediction 
with the enumeration results for the area in the lattice 
model. As can be seen from Fig. \5\ there is good agree- 
ment for small values of J but the data starts to deviate 
away from the predicted curve as J increases. 

When p > Pc{J), the ring is in an inflated state, with 
the area (A) ~ iV^. To obtain an accurate description 
of this regime, we would need to keep higher order terms 
such as R^ /N^ and so on. One thus expects that the lat- 
tice and the discrete problems should differ considerably 
in this regime. 

We now derive expressions for a{J) and /3(J) in both 
the discrete and lattice cases. This is done by considering 
a semifiexible chain subjected to an external force. We 
obtain a perturbative solution for the partition function 
in the limit of small forces. From the partition function, 
we obtain the free energy of the ring. By comparing this 
with the form of Eq. p^ . the values of a{J) and /3( J) 
can be obtained. 



A. Discrete Case 

Consider a semifiexible chain of N monomers. When 
the chain is pulled by a force /, the partition function is 
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given by 



B. Lattice Case 



Z(J,.f.N) 



N 

n 



(18) 



We work in the limit of small forces, treating the J term 
exactly. We consider the / term as a perturbation on the 
zeroth order partition function [ / = in Eq. (|18p ] , given 

by 



Zo{J,N) = [27r/o(J)]^\ 



(19) 



where Io{J) is the modified Bessel function of the first 
kind of order 0. We then expand exp 



t ) as a 



series in / and average each term with respect to the 
zeroth order Hamiltonian. On computing the averages, 
the partition function is obtained as 

lnZ{JJ,N)=liiZo + Nb2f+Nbif + Oif), (20) 

where the coefficients 62 and 64 are given by 



62 = 
64 = 



4(/o-/i)' 

2/2 lo + 3/1 



h-h 



(21) 
(22) 



The InS are modified Bessel functions of the first kind. 
Their J dependence has been suppressed in the equation 
above. 

The mean end-to-end distance in the limit of small 
force is obtained from R ^ d\nZ/df: 



^=2b2.f + 4b4f + 0{f). 



(23) 



For a lattice polygon, where each individual step can 
point only in four directions, we solve the problem of a 
semiflexible chain subject to an external force using the 
exact 4x4 transfer matrix. The transfer matrix in this 
case is given by 



/gJ+Z e//2 e- 



=//2 



gj g-//2 g-J 

e-//2 e-'-f e-f/^ 



//2 



-//2 gJ 



(28) 



We determine the largest eigenvalue up to order and 
hence calculate the partition function: 



lnZ(J,/,7V) = N 



ln(2 -f e-' + e') + —f 



.(29) 



We then follow the same procedure as for the discrete 
case, finding R/N in terms of /, inverting this equation 
to find /, and finally using this expression to compute 
the free energy. We thus obtain 



le-^'^(3e^'^-l) 



i?4 

iV3 



The expressions for a(J) and /3( J) are then 
1 



a{J) 
(3{J) 



in 
1 

— e 
12 



-3J/ 



iSe'-' - 1). 



(31) 
(32) 



MEAN FIELD THEORY 



Solving for / from Eq. we obtain 



f=—— -^f— 
■' ~ 2b2N 46^ [^N 



The Flory free energy F{R) 
to 




(24) 

hiZ + fR, then reduces 



-InZo + 



J_R^ 



bj R* 
1661 7V3 



pR^ 



(25) 



Comparing with Eq. (fT2l) . the factors ck(J) and f3{J) are 
obtained as 



«(J) = 



1 lo — h J^oc 1 



47r /o + Ii 



IBttJ' 



(26) 



/3(J) = 47r2a(J)2 



h + 3Ji 
/0-/1 



2/2 



' 64J' 



(27) 



In this section we present mean-field theories to cal- 
culate the dependence of area on pressure and bending 
rigidity. In Sec. lV Al we address the ISIR model (J = 0). 
The mean field theory presented in [20| performs poorly 
with respect to the Monte Carlo data when p > pc- Here, 
we present an improved variational mean field which re- 
produces the behaviour of the area above the transition 
very accurately. It also yields the correct asymptotic be- 
haviour for the area in the limit of high pressures. In this 
approach, the constraint of fixed link length in treated 
exactly while the closure constraint is satisfied in a mean 
field sense. However, such a mean-field theory fails to de- 
scribe the collapsed phase, also yielding incorrect results 
for the case of nonzero J. 

In Sec. IV B[ wc generalise an earlier mean field theory 
for the freely jointed chain to include semi-flexibility, im- 
posing the constraint of flxed bond length via a Lagrange 
multiplier [Tsl . [20| . The closure condition is imposed ex- 
actly. We thus derive expressions for the average area of 
the ring for all pressures and bending rigidity. 
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FIG. 6: Comparison of Monte Carlo data with the two mean- 
fields for the flexible ( J = 0) case. The density matrix-based 
mean-field approach provides an accurate description of the 
area for p > pc. 



A. Density matrix mean-field for flexible polymers 

In variational theory, a trial density matrix p is cho- 
sen to approximate the actual density matrix[23]. The 
variational parameters are determined by minimising the 
variational free energy Fp with respect to the parameters. 
The simplest mean-field theories assume a trial density 
matrix that is a product of independent single particle 
matrices, i.e. 



P - 



3 



(33) 



where pj is the single particle density matrix of particle 
j. The variational mean- field free energy is 



Fp = {n)p + TY,Trp,\np,. 



(34) 



The variational form for the density matrix should satisfy 
the constraint Tr pj = 1. 

We choose the single particle density matrix based on 
the high pressure limit. In this limit, the ground state 
of our Hamiltonian is a regular N-gon, where the angle 
of the tangent vector is 9j = 27rj7iV. The single 
particle density matrix has a delta function peak at this 
value. At intermediate pressures, we therefore take the 
form of the density matrix to be a gaussian of width a 
(the variational parameter) centered about 9j: 



1 



271 a erf [7r/\/2cr 



exp 



N I 



2a2 



(35) 



where the normalisation ensures that Tr = 1 and erf (a;) 
is the error function defined as 



erf(a;) = —— I e 
A Jo 



dt. 



(36) 
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FIG. 7: Comparison of Monte Carlo data with the two mean- 
field approaches for the case J = 1. 



Using this form of the density matrix, we obtain 



1 , 0Fexp(^V(2a2)) 



2 V2CTerf[7r/V2 



In \/27rcrerf[ 



\/2cr 



,(37) 



where 



J.J^^ _ [(^ " ^^')/V2a] + erf [(^ + ia^)/V2a] 

^ ' 2erf[7r/y2a]e'^V2 ' ^ > 

When 1, the pressure and bending terms in Eq. ([57]) 
can be combined, and the problem is equivalent to one 
of a flexible polymer {J — 0) with an effective pressure 
Pcff = P + J- 

The variational parameter a is chosen to be the a* that 
minimises Fp in Eq. (j37p . This is done numerically. The 
average area, equal to —dFp/dp, is then given by 



{A) = - cot 



47r 



K^a*). (39) 



We now derive the asymptotic behaviour of area in the 
limit of high pressures. We work in the limit when N is 
large. For large pressures, we expect that a* tends to 
zero. In this limit 



K{a) ^ e~'^''/^, cr->0. 
and the variational free energy is then given by 



Fp{a) = N 



(40) 



(41) 



where p — Np/{4tt). Solving dFp/da* = 0, it is straight- 
forward to obtain 



1 



1 - 2 J 



2p 4\/2p3/2 ' 



p oo. 



(42) 
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The area then reduces to 

(A) 1 4J- 1 ^ 

For flexible polymers ( J — 0), this mean-field theory 
reproduces the p > Pc behaviour very accurately. It also 
obtains the correct asymptotic behaviour. In Fig. [6l we 
compare the Monte Carlo data for J = with the results 
of the above mean field theory and contrast it with the 
meanfield theory of Ref . [l^l . 

The density matrix mean-field however, fails to cor- 
rectly obtain the behaviour for non-zero values of the 
bending rigidity. It predicts a first order transition for 
J > 1, in disagreement with results from scaling theory. 
We compare the results of this mean field with the Monte 
Carlo data in Fig.[7]for a system with J = 1. This mean- 
field approach then predicts a transition at p = 0. The 
discrepancy between the two curves increases for larger 
values of J. 

We now describe an alternative mean-field approach to 
this problem which extends the harmonic spring-based 
mean field theory of Ref. [l^l to non-zero values of J. 



B. Harmonic spring mean-field for semiflexible 
polymers 

We follow the approach of Refs. [3, 113] wherein the 
rigid links between particles are replaced by extensible 
springs. The spring constant A of the springs is identi- 
fied with a Lagrange multiplier, chosen so that the mean 
length of a spring equals unity. 

Consider a partition function for N particles given by. 



Z = l dtj exp 



k<j 3 3 



(44) 

Note that while pressure couples to t, the bending rigidity 
couples to the unit vectors t. We make the approximation 
of replacing thyt. This makes the problem analytically 
tractable. 

Expanding the tangent vectors in Fourier space as. 



Y^[Ak cos{jk) + Bk sin(jfc)]. 



\[^Y}^k cos(j/c) + B; sin(jfc)], (45) 



By completing the squares, this integral can be written as 
a gaussian integral and hence can be calculated exactly. 
This gives 



n 



1 



A — J cos k 



P 



4fc2(A - Jcosky 



(47) 



The parameter A* is determined by equating the mean 
square link length to one, i.e 



1 dlnZ 
N dX 



1. 



(48) 



This gives 



N = 



N 

E 

1=1 



A* 



Jcos(M) 



2f 



P[X*-Jcos{^)]^-p 



(49) 

where p — pN/An. 

When J — 0, the first factor in Eq. ((49)) becomes in- 
dependent of /, and then the resultant expression can be 
evaluated exactly. Hence, an analytic expression for A* 
can be obtained in this case [l^l. For J ^ 0, this is no 
longer possible, and for finite system sizes the resultant 
equation must be solved numerically. When N ^ 1, it 
is still possible to extract the behaviour of the system 
analytically. 

We now determine the phase boundary from Eq. (I49p . 
We will consider the limit N ^ 1. First, note that A* — 
Jcos{2Trl/N) for all I. For positive A*, this gives the 
condition that A* > J. Second, consider the term in the 
denominator for ^ = 1. It is (A* — J)^ —p^- If we assume 
that A* is continuous in p, we have the second constraint 
that \* > J + p. 

Setting X = and converting the first sum in Eq. (I49[) 
to an integral, the equation for A* reduces to 



1 = 



1 



1 



VA*2 - J2 iV(A* - J) 



2k 



The sum in Eq. (|50p is convergent if the ratio p/ (A* — J) < 
1. In this case, we keep only the first term on the right 
hand side of Eq. ((50)) . This gives. 



A* = Vl + J^, for p < 



(51) 



The critical pressure is obtained when the ratio p/ (A* — 
J) becomes equal to 1, i.e. 



where k = iirl/N , I — 1, 2, • • • ,N. The partition func- 
tion then reduces to 



^ = n / dAkdA^dBkdB'^ 



(46) 



,( J) ^\* -J = vT+J^ _ J. 



(52) 



For large values of J, this goes as Pc{J) ^ 1/2J, which 
differs by a factor of 2 from the answer obtained by scal- 
ing arguments [see Eq. ((H))]. 

We shall now estimate A* in the different scaling 
regimes. We assume that A* is a non-decreasing function 
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of p (as in J = 0) . Then, since we have the constraint of 
A* > p + J, the ratio p/ (A* — J) must continue to remain 
at 1 for p > pc- Thus, above the critical point, we obtain 


3 


J = 0.5 ^^^^ 


X*~p + J, for p > Pc- (53) 

However, a simple substitution of Eq. (|53p in Eq. (|49p 
for p > Pc does not satisfy Eq. . We therefore need to 
calculate the correction term arising from large but finite 
A''. We start by considering Eq. ((49l). The first term can 
be summed exactly, giving 


^ 2 

s 

1 





1 - 



VA2 



N 



J2 
1 



2f 



Ea- Jcos(2^)P(A- Jcos(^)) 



p^ 



.(54) 



We calculate the finite size corrections to A* as follows. 
Let 



A* 

P>Pa 



p + J-6. 



(55) 



When (5 — > 0, the main contribution to the left hand side 
of Eq. ([5^ comes from the I = 1 term. The contribution 
from other I is convergent as (5 — > 0. Expanding the right 
hand side as a series in 6, we obtain 



N 



1 



Sip + J) 



+ 2pJ (p2+2pJ)3/2 



(56) 



The S independent term in the right hand side of Eq. (|56p 
is nonzero for p > pc and is equal to zero for p = pc- 
Thus, when p > pc, we keep only the first term in the 
right side, while at p = Pc, we need to keep the second 
term too. Solving for S, we obtain. 



1 

yjV ( 1 + J2)1/T . 
^p2+2pj-l' 



N 



P^Pc 
P> Pc 



(57) 



We are now in a position to calculate the mean area 
(A) from This gives. 



(A) 



Np 



N 



1 



2tt ^ - Jcos(- 



N 



(58) 



The numerical values obtained for A are then substituted 
in this equation to get the corresponding value of the 
area. We can, however, analytically determine the scaling 
behaviour of the area in the limit of large system sizes 
from the values of A calculated above. 
For p < Pc, we have. 



N 



1 



P{y'l + ,P- Jcos{2ttI/N)Y - p'^ 



(59) 



At the critical point, we obtain, from Eqns. (j56p and 

p^Pc- (60) 



Air 



-0.5 







0.5 



FIG. 8: Area collapse for flexible and semiflexible polymers 
around the critical point. This verifies Eq. (|62|l . The data is 
for iV = 80, 100, 120, 140, 150 for the lattice problem. 



Similarly, for pressures greater than the critical pressure, 
we obtain, from Eqns. ([55)1 and 



A^ 
47r 



1 - 



1 



Vp^ + 2P J 



->oo 1 

p 



J 

2^ 



P> Pc 



(61) 

This mean field theory reproduces the qualitative be- 
haviour of the simulation data correctly. It predicts a 
continuous transition for all J, unlike the density matrix 
field theory. However, there is a quantitative disagree- 
ment with the data. This can be seen by comparing the 
results of this mean-field theory with the simulation data 
in both the flexible (Fig. [6]) and semi-flexible (Fig. [7]) 
polymer cases. 



VI. 



SCALING AND CRITICAL EXPONENTS 



The order parameter that describes the collapsed to 
inflated phase transition is the ratio of the area to the 
maximum area. When A^ ^ 1, the ratio is zero below 
the transition and non-zero above it. The behaviour near 
the transition line can be described by the scaling form 



(A) 
A 



(62) 



where 0, (3 are exponents and g{x) is a scaling function. 
When X — > 0, then g{x) — > constant. When x oo, 
then g{x) ~ x^ . When a; ^ — oo, then g{x) ~ l/x [see 
Eqs. and p6|) ]. This immediately implies that 



(1+/3) = 1. 



(63) 



To obtain the one independent exponent, we resort to 
the scaling theory (see Sec. ITV)) . At Pc, (A)/Amax ^ 
1/\/N. At the critical point, the area scales as N^^. 
Combing with Eq. ([55)1 . we obtain 0=1/2 and (3=1. 
These exponents are independent of J. 
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FIG. 9: The scaling plots for compressibility x when scaled as 
in Eq. (|66|l . The data is for the lattice model with J — 0.5 and 
J = (inset). The system sizes are = 80, 100, 120, 140, 150. 



FIG. 10: Collapse of the En{A) for different values of A'^ when 
plotted against /N^. The data is for the lattice model with 
J = 0. 



In Fig. [51 we show scaling plots when area is scaled as 
in Eq. ((62| with (j> and /3 as above for the cases J = 
and J — 0.5. The excellent collapse shows that the Flory 
type scaling theory gives the correct exponents. 

We now look at the fluctuations. Consider the com- 
pressibility X defined as 



X 



An 



1 d{A) 
dp 



(64) 



When p < Pc, X can be calculated from Eqs. and 
(flB to be 



X 



1 



p^ p2sin^(7rj5/pc) 



, P<Pc 



(65) 



Thus, X diverges as (pc—p)^^ below the transition point. 
The behaviour near the transition point is described by 
the scaling form 



X~ Af^^/i[(p-Pc)iV^] 



(66) 



where h{x) is a scaling function and = 1/2. When 
X ^ 0, then h{x) constant. When |a;| ^ 1, then 
h(x) ^ x^"' . Comparison with Eq. (|65p gives 7 = 2. 

In Fig. [HI we plot the compressibility scaled as in 
Eq. (|66p for two different values of J. A good collapse 
is obtained again showing that the Flory type scaling 
theory gives the correct exponents. Similar, but noisier 
data can be obtained for the discrete model. We thus 
conclude that the introduction of semiflexibility does not 
affect any of the exponents describing the transition. 



VII. 



THE LATTICE PROBLEM 



In this section, we present some additional enumera- 
tion results for the lattice problem. Consider the scaling 
theory presented in Sec. IIVI The inextensibility of the 
polymer was captured by the R^/N^ term for a polymer 
of extent R. This was obtained from a calculation based 
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FIG. 11: The asymptotic behaviour of area in the limit of 
large p as computed for the lattice model. The curves are 
straight lines when plotted against 



on the extension of a polymer under a force. Here we 
present numerical evidence supporting this. 

Let P/v(j4) be the probability (at P = 0) that a walk 
of length N encloses an area A. In Appendix El we show 
that [see Eq. (|M1) ] 



Pn{A) = 1/ , ^, TV ^ oo, ^ fixed. 



where the scaling function I{x) is given by 



I{x) ~ TT sech^ (27ra;) . 



(67) 



(68) 



We consider the corrections to the scaling form in 
Eq. dnH). Let 



En (A) 



NPNjA) 
I{A/N) ■ 



(69) 



Scaling theory predicts that En (A) should be a function 
of one variable A^/N^. This is verified in Fig. [TOl where 
InEjyiA) is plotted against /N^ for a range of system 
sizes. 
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We also study the behaviour of area when p is very 
large. When j5 ^ 1, the behaviour is seen to differ from 
the discrete version of the problem. It can be shown to 
beH 



1- 



A 



1 

— , P 



(70) 



This should be contrasted with the discrete case which 
varied as 1/p. In Fig. 1111 we show numerical confirmation 
of the prediction of Eq. ((70)) . 



VIII. 



CONCLUSIONS 



In this paper, we have proposed and studied lat- 
tice and discrete models for self-intersecting pressurised 
semi-flexible polymers. Our work generalises results of 
Ref. 20] to include a bending rigidity. A simple vari- 
ational mean-field approach provides very accurate fits 
to the Monte Carlo data for this problem in the absence 
of semi-flexibility. The mean-field approach for J — 
was generalised to the semiflexible case. The 
phase boundary between collapsed and inflated phases as 
well as expressions for the area as a function of p and J 
in the different phases were obtained analytically. 

We have shown that the essence of the physics is cap- 
tured through simple Flory approximations. The scaling 
predictions of the Flory theory were verified numerically 
for both the lattice and discrete cases. 

We have also investigated the behaviour of the system 
in the extreme limits of a fully pressurised polymer ring 
and a collapsed configuration. For the fully pressurised 
ring, we deduce the leading order asymptotic behaviour 
of the area in both the discrete and lattice cases. The col- 
lapsed phase was studied by mapping this problem onto 
a quantum mechanical problem of an electron confined to 
two dimensions and placed in a transverse magnetic field. 
The analytic results thus obtained fit the data accurately. 

The usefulness of these results for more realistic sys- 
tems lies in the fact that both the restriction to the signed 
area as well as allowing for self-intersections at no en- 
ergy cost are irrelevant in the large p limit. The results 
obtained at large p should therefore apply both quali- 
tatively and quantitatively to the more realistic case of 
a pressurised self-avoiding polymer, where the pressure 
term couples to the true physical area and not to the 
signed area. This is the LSF model Q. The approach 
presented here is thus also useful in imderstanding the 
behaviour of a larger class of models, some of which are 
more physical in character, but which lack the analytic 
tractability of the model proposed and studied here. 
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APPENDIX A: ANALYTIC ANSWER IN THE 
SMALL PRESSURE REGIME 

The problem of self-intersecting polymers in two di- 
mensions with no bending rigidity {J = 0) is analogous 
to the quantum mechanical problem of an electron mov- 
ing in a magnetic field applied transverse to the plane 
of motion. Using this analogy, we obtain analytic ex- 
pressions for the partition function Z and C'n{A), the 
number of closed walks of area A. 

When an electron goes around a magnetic field, it picks 
up a phase factor proportional to the flux enclosed by 
the path. This flux is proportional to the product of the 
strength of the magnetic field times the algebraic area en- 
closed by the loop. The propagator then is the sum over 
all such loops. This suggests that the partition func- 
tion for the polymer problem can be obtained from the 
quantum mechanical propagator for the electron problem 
provided the constants are appropriately mapped. 

For an electron of charge e and mass m in a constant 
external magnetic field B, in the z direction, in the case 
when the electron returns to the origin, the kernel can be 
written as 12511. 



^ ' ' ' ^ \2mhtJ Vsinwt/2 



where m — eB/mc. It picks up a flux $ given by 

eBA 



he 



(Al) 



(A2) 



The motion of a quantum mechanical particle is gov- 
erned by the two-dimensional Schrocdinger equation, 

The classical diffusion equation for a particle in two- 
dimensions is 



dP 1 fd'^P d'^P 
4 \ dx^ dy'^ 



dt 



(A4) 



where y) is the probability of finding the particle at 
{x,y). Thus to map the results of the quantum problem 
onto the polymer problem, we need to map t —it and 
also identify 



A - 1 
2m ~ 4' 

Also, to match the diffusion constant, we need, 
leB 



he 



(A5) 



(A6) 



Substituting in the propagator and equating it to the 
partition function, we obtain 



P 



(A7) 
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where we have introduced the appropriate normaUsation 
factors. C'pf{A) is now obtained from the partition func- 
tion by performing the inverse Laplace Transform with 
respect to p. This gives 



Cn{A) 



■sech 



,27rA 



and 



^ = r cot -r 



(A8) 



(A9) 



The free energy wih have a singularity at p = Att/N . 



Below this p, the expressions are valid for both the dis- 
crete case and the lattice. Exactly the same expression 
has been obtained by using the harmonic spring approx- 
imation 17]. The expression for area matches both the 
simulation and lattice data quite closely for low pressures, 
as can be seen from Fig. H) 

Moreover, if we recall the Flory prediction that by 
rescaling area and pressure by Pc{J), we can obtain the 
results for non-zero values of the bending rigidity from 
the answer of the problem with J = 0, we see that 
the above analysis also predicts the area expression for 
nonzero values of J. 
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